home *** CD-ROM | disk | FTP | other *** search
/ Libris Britannia 4 / science library(b).zip / science library(b) / ASTRNOMY / AA_51.ZIP / LIGHTT.C < prev    next >
C/C++ Source or Header  |  1993-02-09  |  2KB  |  79 lines
  1. /* Correction for light time from object to earth
  2.  * including gravitational retardation due to the Sun.
  3.  * AA page B36.
  4.  */
  5.  
  6. #include "kep.h"
  7.  
  8. int lightt( elemnt, q, e )
  9. double e[], q[];    /* rectangular position vectors */
  10. struct orbit *elemnt;    /* orbital elements of object q */
  11. {
  12. double p[3], p0[3], ptemp[3];
  13. double P, Q, E, t, x, y;
  14. int i, k;
  15.  
  16.  
  17. /* save initial q-e vector for display */
  18. for( i=0; i<3; i++ )
  19.     {
  20.     p0[i] = q[i] - e[i];
  21.     }
  22.  
  23. E = 0.0;
  24. for( i=0; i<3; i++ )
  25.     E += e[i]*e[i];
  26. E = sqrt(E);
  27.  
  28. for( k=0; k<2; k++ )
  29.     {
  30.     P = 0.0;
  31.     Q = 0.0;
  32.     for( i=0; i<3; i++ )
  33.         {
  34.         y = q[i];
  35.         x = y - e[i];
  36.         p[i] = x;
  37.         Q += y * y;
  38.         P += x * x;
  39.         }
  40.     P = sqrt(P);
  41.     Q = sqrt(Q);
  42. /* Note the following blows up if object equals sun. */
  43.     t = (P + 1.97e-8 * log( (E+P+Q)/(E-P+Q) ) )/173.1446327;
  44.     kepler( TDT-t, elemnt, q, ptemp );
  45.     }
  46.  
  47. if( prtflg )
  48.     printf( "light time %.4fm,  ", 1440.0*t );
  49.  
  50. /* Final object-earth vector and the amount by which it changed.
  51.  */
  52. for( i=0; i<3; i++ )
  53.     {
  54.     x = q[i] - e[i];
  55.     p[i] = x;
  56.     dp[i] = x - p0[i];
  57.     }
  58. showcor( "aberration", p0, dp );
  59.  
  60. /* Calculate dRA/dt and dDec/dt.
  61.  * The desired correction of apparent coordinates is relative
  62.  * to the equinox of date, but the coordinates here are
  63.  * for J2000.  This introduces a slight error.
  64.  *
  65.  * Estimate object-earth vector t days ago.  We have
  66.  * p(?) = q(J-t) - e(J), and must adjust to
  67.  * p(J-t)  =  q(J-t) - e(J-t)  =  q(J-t) - (e(J) - Vearth * t)
  68.  *         =  p(?) + Vearth * t.
  69.  */
  70. velearth(TDT);
  71. for( i=0; i<3; i++ )
  72.     p[i] += vearth[i]*t;
  73.  
  74. deltap( p, p0, &dradt, &ddecdt );  /* see dms.c */
  75. dradt /= t;
  76. ddecdt /= t;
  77. return(0);
  78. }
  79.